


Institutional Archive of the Naval Postgraduate School 





Calhoun: The NPS Institutional Archive 
DSpace Repository 


Theses and Dissertations l. Thesis and Dissertation Collection, all items 


1992-12 


Computational and experimental investigation 
of the aerodynamic characteristics of a 
windsurfing sail section 


Avila, Matthew R. 


Monterey, California. Naval Postgraduate School 
http://hdl.handle.net/10945/23577 


This publication is a work of the U.S. Government as defined in Title 17, United 
States Code, Section 101. Copyright protection is not available for this work in the 
United States. 


Downloaded from NPS Archive: Calhoun 


Calhoun is the Naval Postgraduate School's public access digital repository for 


۹ DUDLEY research materials and institutional publications created by the NPS community. 
FW ο ος Calhoun is named for Professor of Mathematics Guy K. Calhoun, NPS's first 


WW KNOX appointed — and published — scholarly author. 

“Ἢ LIBRARY Dudley Knox Library / Naval Postgraduate School 

411 Dyer Road / 1 University Circle 
Monterey, California USA 93943 





http: //wwwenps.edu/library 





[ ایا لها اد ۱ه 
Mw rima 1 ut ۱‏ ` 
n 1 AUN ey + TT Me RX | ۲ YA y 24)‏ 17 
ia rtt 0‏ رق E : v ۱ | ۱ ' 7 À‏ 
«T Y M or λα... τη‏ 
Dn E E A uh‏ 
j EX!‏ نام : | ۱ ۱ Y ας d 1 1 Di,‏ ) 
P 4 a:‏ 
N i ۱ ۱ Y i ۱ a ' ur u D‏ را 5 ή‏ 
i T 1 u ON M NEC Y | ' Γι : M‏ 1 1 
er b ` L ΙΙ. T LN Ts‏ 
i PC‏ یب ا اد U Eu‏ : 
y ONES eh, “iv xa‏ 
r . à ۴ 1 4 1 F $u 8 " " =“‏ 
ΠΡΙ, MNT aur n ος ١‏ اام να os pa‏ 1 
" اا لس | DA. T AE a mb TP u i‏ 
l Lè B ET ۱ Fd (n y 1‏ 
A AP E f j ۸۰‏ یڈ A 7 i‏ ۰ 
À ٣ y IP / [ 1 IA ΠΤΙ y 1 AE A‏ : 
ῃ k i 1 I 1 1 i 10 [ 1. 41} 31‏ 
E ٩ a σος N TA‏ 
EE μη 1 i ui cM‏ 
ا e‏ 1 زا LI 1 4 ΠΠ ur‏ 
P ! e Pai ' 1 Ara Wow "P" yt‏ 
i A MAA 17 1 :‏ 1 ۷ ال j‏ 
n L 7 I y Pw 8 fi % y : 0‏ 
a 1 ¡de i ۹ M κ M ee a |‏ لك 1 
t F η P out 1 "TR X 0‏ ' 
u‏ پل ῃ IM^ i ٨ Na ۱ 4 +! ۷ RA 4 NS A Ns‏ 9 
σε S n^; i QUA μην EUN‏ 
D ? * dg E * A | la k ^ 3‏ : : 
P * TA D M n τς E 1 Kv = η‏ " اي ۹ N E N n I DOLI ET ^ ۱ 1 SI‏ 
ἃ E T OT a li : e d‏ رب d ٩ L 7 2 ` Ll‏ 
مہا رام ال می T. M‏ وه ې ٦ à‏ - ! 1 
$ ړا Y ee K M A‏ تی یں ات ae 7 | f y ' A m A I a . 1 e M‏ ۱ : 
D 4 [NE k 4 n L 5 4 [ My A Y X» T bie 4 0‏ 5 
L , |: is | e HL ۱ 1 4 ¿Ne 9 aa Ag W iv ae y Í y xu 1 ΤῚΣ Tama MS DET S oa ἘΣ 4‏ 
hi‏ پد 1 5 0 A ' ή‏ 
αν E uo ow, AUN. 7 M To ۴ y m s N X τ - py DS Mu RER‏ ۱ 
P‏ کک y D D " ۰ LEO 7 À F | "T! % : xl m [ * ii 4 Ἢ ur pv) CCS‏ ‘ 
AA ek: η‏ کت i ae j N^ E T ru P M MN. d E [ 4 x ۷ C 17 E kna‏ : 
٢ | r ۸ 5 and 108 : Ἂ A $t AN‏ و Ç y f M : , A Uu‏ | 
D D " 1 : : ] 1% ^ "m e 7 7 ή Er‏ 
AN ۳ KE‏ پا I Atte as ME ü n3 "YT I A e ο το αὐ, i Py DOM‏ 
iN ib. cn | Ww E d FUA‏ ا0 (i! ٩‏ 1 1 1 - 
o Y‏ ا لے کےا سس y‏ و ۱ : پا » : : 
U * ٠ : : ab L de ñ ٩ % ^ : I oe e 3 A‏ و , u y i‏ 3 : .: 
Bv. H vale > I τ ks EI N ΚΝ Ny CH n T‏ در ῃ d A‏ 
ا ۹ Í‏ پک ۱ 1 Li . σος "mE X. ` 7 rf ١ h 354 " J ARM F)‏ و" A L‏ 
i Al |‏ و چا a E 4 J y M‏ 7 2 د : | Je 0 y‏ 
ποι. : 1 ' ۳۰ κ i^ E s Wave B PA AR mx.‏ 
E 5 1 4. ἡ Dels UE * n s AIDA i‏ . : 
E‏ اځ < A‏ .70ء : ] ῃ‏ 
SSe, 8 I AM Ns Ww 7 4 SA dE i pe I κ ۱ : : ha‏ ا یہ 1 ` : E i 1 1 f‏ " 
XO s go, 3 7 1 " A "WW Y igh ' A CN LO R 7 tr s. ۹ WA iN κ) δον‏ 
a, ` i Lj i L 1 7 E k * ΄‏ 1 
Ἵ N 4 X B i μα I `. - <‏ 
m A Ἴ n A i a da RA ^w | , 4 " " ٧ N ER‏ 2 
y 1 » 22 2 ı$2 E 4 y " I tu "a CN ar‏ 0 
ἐπ i ۱ αι‏ 1 وس 1 M ΠΠ; ΠΠ‏ 7 
tp KANT αν oer Me V. ERA Y‏ , 
په : ہل A A LAMA‏ ء7۹ a‏ 
y a a‏ ا j 3 ' tud A n T‏ 
A 1 N 0 r P ο 0 7٦ W s πὶ‏ - ۰ 
i CE, «git E. d EC | 1 V‏ ; : 
5 ا 0-0 ; یہہ e τ m JL δ | Ti ٠ | a a m‏ 7 ] ` 
t ae, 1 5 re IRS Το = 1 E‏ 20 و " f 4 : 1 νο‏ زا سه | : i D 1 x y r‏ 
XE AA IA Π uA κῶν τν μὴ Punt LO‏ : ͵ 0 : 1 ' 
E 4 J 1 A I 1 1.4 E y a οἷ, ος ۳ po w η‏ 
"m H 1 j D‏ " 
RR i ٢ M Ν 0‏ یله چ 27 S 8 E‏ 
T": ΓΣ ΠῚ 1 A ^ 7 bey 43 NOM P E Í‏ 
o 1 [ 5 1 Pal 3 š‏ 
S es‏ رر ἢ D a I A = * 2 ۶ h^ πον.‏ 
a " E] 7 D W ٢ ۳ Y [ νων‏ 
ل۷ ۱ 0 i | : d i 1 i 0 Ln : i Lar 0) D‏ 
r tow E a‏ 1 : " 
M 1٩ i 1 E |‏ " 1 ? ` 
E n 1 d / I y " P en‏ 0 
)و u‏ اع ٤ : i‏ 
ἡ, n m CA 5 CA u‏ له Y MEE‏ , 
F 7 ^a Β š‏ 
y : a i d 1 A : 1 i‏ 1 
A J 4 j 4 i 7 ya 8 ۱ eo a‏ 
euer ; ^ 1 A "P: Aj ΠΠ 0 P ۱ M‏ 
z i I ε8 P " 7 1 : L‏ ام 1 را Ν : = i‏ ہیں i‏ 
٢‏ صا 2 1 A UR 1 x‏ 
«a i k ٦ i ^ " ۹ : ۱ ' vs‏ » 
k 9 V 7 A E. A ! P $ 1 κ‏ ۱ 
u ! 4 ἘΣΤΙ ΕΚ‏ 1 
D = H 1 1 / 3 i i * | é‏ 
he ^‏ و ۳ t‏ " : 
y , |‏ ' 
ZW‏ : ہی a | "e ١ ۸ Oy‏ ۰ .. 
A νο : 7 i ΠΣ $ ` ۷ E‏ 
i ` I ' D M * = TE‏ 
y 1 „ui 4‏ 
à I i -‏ : ` 
Ü : 4 7 A Π 1 j LI 1.8 1 i í Th 7‏ 
I 4 I‏ 1 
LÀ - i 4‏ ^ 1 دا Au > i = A‏ 1 ۰ 
V n - LI kato ` * EN‏ ۰ 
I 7 Š kE O ὁ Ms B Y m‏ " 
MO E I | aM / ۱ 9‏ | 
DN ۲ : A b P. p 5 7‏ * " : ۰ 1 
d 4 [ σα N 2‏ » 
i y y‏ ) | .. ` " 
١‏ 7 
I κ ; 1 B‏ 
M ۱ i i |‏ : ' 
j 1 Li η .-'‏ 
P è : | : e I S ۱ ad a‏ 
η i $ 2 i 34 | r E P‏ 
j ۰ = j 1 7 D š 5} ir‏ 
E i ۱ : 1 , . N y a i‏ 
g i J | N : [| 1 y‏ : 
E ۱ 1 -‏ 
ra |‏ 1 
E M E M d‏ 3 , 
i K fi L I E L4 j M ۱‏ 5 3 
T NET - : A 7 E e‏ $ 5 ۱ 
E " 1 i ' m‏ 
η bs : 1 1 M : =‏ 
ا نے 8 i s ` 7 / ! LA 1 3 a‏ 
Fo ! - 2 À‏ 3 
2 1 ' 
ῃ , n ` 1 E Vm,‏ . 
I 7 ١ Sr] κ " ri i ,‏ 
n : I ' M P ` E mi ' i ` Y i‏ ' 7 
Ἢ k dı N 1 y yi -‏ 7 
سا 1 : D‏ 
E k, 7 (> NE an‏ 
D , E ` ñ L 7 [| / ^ -qm‏ 
S Di Re‏ ۲ھ : : - 2 Li‏ 
a Koa 3, = . =F‏ " - ' 
I I =‏ - 
i r L : s‏ 
٤ م١ ë [ n da "T Td‏ + 
I ο E | a d yt‏ 
n - l d I pi | " i 71,28]‏ ; 
i‏ 
U ۴ -‏ 
O I : "-— i‏ 0 
ای 17 7 
داو للا TP : Ξ : I : ΛΙ‏ 
q " " [ EM si‏ 
i >‏ | 4 
| 4 ' ۰ 
δ 1 ٩ i 1 i 3 » ¡e κα !‏ 
A d à‏ 
LI : >‏ 
i ١ y‏ 
Y 4 m Dd TEE BETA‏ ^ 
i D k -‏ 
a " 1 EA C 1 H‏ 
ç Ye š: MW : M Lift‏ 
N ή (y : -i Y ΙΙ. tei ) Se‏ 
٦ ar - ¡E‏ 4 
E το d‏ 5 کر 6 d BE n‏ " و 1 
يب ΓΚ LE‏ 0 [ 1 - = 
E] l 1 i - [| : one‏ 0 $ 
if‏ " 7 " 
M N ^. va 4 1 1 1 SE‏ 
h | E 1 4 Ἡ * ٩ E‏ 
H m ۱ I 1.3 5 1‏ 
P on f "m‏ کی رر i‏ , 
A‏ 11 اا کا | 4 
š M hi > 4 : ۷ n "aT 8‏ 1 
iy ٠ ΜΙ Ç J Hu Y a i‏ ۱ 
A $ š Ἢ f eur PEE ON‏ 
r] 4 1 D" L 1 nh ^ PE" ° 1‏ 
I u4 Ca d Fi I ral a‏ ۱ 0 
E i E ON AME IN ΄ απ 1 E. 7I‏ 
یپ i 1 si R K 3 E]‏ 
E [ 8 :‏ | 
ay Boxer wurde T‏ 
h E - + À UN.‏ 
EU‏ ا " و M an‏ 
i is ٨ en ean αλ ΝΟ.‏ 
0 " « | > 
I i LA ۱ / M‏ 
i i , mo As - Ki /‏ ۶ 
“N [2 y 5 : 1 L ot |‏ 1 
E “u. : - - s‏ 1 
A " 1 *‏ دی ; y‏ سے LE‏ 1 
ary η ΓΝ‏ : | 
ma .‏ اف T +? i ۱ να ER‏ 1 
کر ^ 1 0 I a‏ : 
j. a ui TARI‏ 
E ۱ E u 1 -‏ يم à‏ 
E E.‏ د = 1 
š | ñ " IN 1 2 i] A ۱ 2 T į ۳‏ 1 - 
αὶ : i V dei Y 7 ; y a 02‏ 1 
I‏ - = ][ دو || i U‏ + " > 
Uu j a | y LIE á dues‏ % لن I‏ 1 ` 14 
١ T xli 1 ۶ C ΣΝ‏ 
|( اند i= i T [E i er‏ 7 
ΚΑ. Πω‏ و e‏ 
h‏ : 0 دې : : I‏ 
ade,‏ سر لیا d 1 l LS AX. αλ ο δ xs‏ : 
4 | " 

















Unclassified 
Security Classification of this page 





REPORT DOCUMENTATION PAGE 


la Report Security Classification: Unclassified 15 Restrictive Markings 

2a Security Classification Authority 3 Distribution/Availability of Report 

2b Declassification/Downgrading Schedule Approved for public release; distribution is unlimited. 
4 Performing Organization Report Number(s) 5 Monitoring Organization Report Number(s) 

64 Name of Performing Organization 6b Office Symbol 7a Name of Monitoring Organization 

Naval Postgraduate School (if applicable) *52 | Naval Postgraduate School 

6c Address (city, state, und ZIP code) 7b Address (city, state, und ZIP code) 

Monterey CA 93943-5000 Monterey CA 93943-5000 


88 Name of Funding/Sponsoring Organizalion 6b Office Symbol 9 Procurement instrument Identification Number 
(if applicable) 


Address (city, state, and ZIP code) 10 Source of Funding Numbers 


Program Element No Work Unit Accession N 


Il Tithe (include security classification) COMPUTATIONAL AND EXPERIMENTALINVESTIGATION OF THE AERODYNAMIC CHARACTERISTICS OF A WINDSURFING SAIL SECTION 


12 Personal Author(s) Matthew Avila 


13a Type of Report 13b Time Covered 14 Date of Report (year, month, duy) 15 Page Count 
Engineer's Thesis From To macember 992 75 


w 


16 Supplementary Notaíion. The views expressed in this thesis are those of the author and do not reflect the official policy or position 
of the Department of Defense or the U.S. Government. 


18 Subject Terms (continue on reverse if necessary and identify by block number) 


Field Aerodynamics, Computational Fluid Dynamics, Sailing 


EMEN e Ον -........., 

inane 

'19 Absiract (continue on reverse if necessary and identify by block number) 

In this thesis results of a computational and experimental investigation of the aerodynamic characteristics of a sail section used in 
windsurfing sails are presented. State-of-the-art computational methods (panel, direct boundary layer, viscous-inviscid interaction, 
Euler, and steady/unsteady Navier-Stokes) were used to predict the aerodynamic loading and stall characteristics. These 
predictions were found to be in satisfactory agreement with tuft and smoke flow visualization experiments carried out in the Naval 
Postgraduate School low speed wind tunnel at a Reynolds number of 800,000. Further, all computational work was completed on 
the Silicon Graphics Indigo workstation to demonstrate that only modest computer facilities will be necessary for these methods to 


migrate to the field of sail design. 


0 Disiribution/A vailability of Abstract 21 Absiract Security Classification 
_ unclassified/unlimiled __ same as report __ DTIC users Unclassified 


22b Telephone (include Area Code) 
656-2058 








22c Office Symbol 
ΑΑ/ΡΙ 


?2a Name o; Responsible Individual 
Professor M. F. Platzer 










D FORM 1473,81 MAR 83 APR edition may be used until exhausted securily classificalion of this pa 


lL other editions are obsolete Unclassifie 


i 1260188 





Approved for public release; distribution is unlimited. 


Computational and Experimental 
Investigation of the Aerodynamic 
Characteristics of a Windsurfing 

Sail Section 


by 


Matthew R. Avila 
Lieutenant, United States Navy 
B.S. United States Naval Academy, 1984 


Submitted in partial fulfillment 
of the degree requirements for the degree of 


AERONAUTICAL AND ASTRONAUTICAL ENGINEER 


from the 


NAVAL POSTGRADUATE SCHOOL 
December 1992 


ABSTRACT 

In this thesis results of a computational and experimental 
investigation of the aerodynamic characteristics of a sail 
section used in windsurfing sails are presented. State-of-the- 
art computational methods (panel, direct boundary layer, viscous- 
inviscid interaction, Euler, and steady/unsteady Navier-Stokes) 
were used to predict the aerodynamic loading and stall 
characteristics. These predictions were found to be in 
satisfactory agreement with tuft and smoke flow visualization 
experiments carried out in the Naval Postgraduate School low 
speed wind tunnel at a Reynolds number of 800,000. Further, all 
computational work was completed on the Silicon Graphics Indigo 
workstation to demonstrate that only modest computer facilities 
will be necessary for these methods to migrate to the field of 


Sail design. 
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I. INTRODUCTION 

The field cf computational fluid dynamics (CFD) has come of 
age over the past several years. Rapid advances in high speed 
digital computers, their subsequent lower cost and the maturation 
of CFD codes suggest that CFD methods should spread to areas beyond 
the aerospace industry in the near future. The discipline of 
competitive sailing is a promising area that could benefit from the 
use of CFD techniques. Some have suggested that such applications 
would be trivial. However, experience shows that sailing is a 
fiercely competitive sport ranging from amateur levels to Olympic, 
professional and the America’s Cup levels. In practical terms, the 
bottom line has always been performance, which directly relates to 
the incorporation of technology into the design process. The aim 
of this research will be to demonstrate that CFD techniques may be 
readily adapted to the field of sail design. 

From an engineering standpoint, sailing craft are challenging 
due to the interaction or coupling of two different fluid media, 
the water’s surface and the atmosphere. Additionally, the two 
media are subject to separate noise or turbulence in the form of 
waves, wind gradient and gusts. The objective of this thesis will 
be to apply CFD techniques to the analysis of the above water 
portion of the problem and to provide some experimental 
verification. The scope of the analysis will be limited to a two- 
dimensional approach. A full three-dimensional analysis would have 
been optimum, but the complexity of the problem and computing power 


required would have become excessive. 


There are numerous rig and sail configurations which ought to 
be analyzed. Planing hull sailing craft, in general, are more 
challenging because they are capable of higher speeds than 
displacement hull yachts. Planing hull craft, moreover, experience 
a pronounced rise in drag prior to transitioning from a 
displacement hull mode to a planing hull. 

For this study a sailboard or windsurfer rig was selected over 
other rig and sail combinations for several reasons. The primary 
reason for the selection of this rig and sail combination is the 
rig's simplicity, Figure 1.1. Its shape is approximate to a simple 
wing or foil shape of moderate aspect ratio and is favorable when 
compared to a sloop or other multi-sail configuration. The 
acceleration through the above mentioned planing transition is of 
high interest for study. This acceleration or quick planing 
ability is of paramount importance to the competitive sailor on any 
planing hull dingy or sailboard. To accelerate quickly the sailor 
in effect must generate large lift coefficients from the sail. 
This is achieved through a large sheeting angle (high angle of 
attack) in combination with the highly cambered shape built into 
the sail section and in the case of the board sailor through the 
dynamic sheeting or pumping of the sail. The effects of pumping 
are well known to all competitive sailors. However, this sailing 
technique is strictly prohibited in all classes except some types 
of boardsailing. The pumping method utilized is a refined art to 
many board sailors to the degree that sail designers and 


competitors note this quality in sails. 





Figure 1.1 WindSurfing Sail 


The second regime in which the sail performance is critical 
occurs after the craft has fully accelerated on a plane. Asa 
craft accelerates the apparent wind shifts forward from the 
sailor's frame of reference while the true wind velocity remains 
constant, Figure 1.2. While sailing with the apparent wind far 
forward the sail is trimmed at or near the maximum lift to drag 
ratio. The maximum lift to drag ratio in fact limits how close to 
the wind a vessel may sail or the minimum apparent wind angle that 
may be sailed. 

The technical challenge for analysis in this regime of 
aerodynamics lies with the low Reynolds number (Re) below one 
million. The bulk of data to date for airfoil sections is for 


Reynolds numbers above one million. CFD work has been done at 
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Figure 1.2 Apparent Wind and L/D Relation 


Reynolds number below one million, but the work has been frequently 
restricted to internal flows in turbines or ducts or to research 
for high altitude aircraft where the Mach numbers are significant. 
The low Reynolds number implies that the viscous forces will be 
large and require thorough investigation. This will require the 
use of viscid codes in addition to simplified inviscid codes to 
understand the nature of the flow field. On the positive side, the 
problem is somewhat simplified in that it is entirely 
incompressible. Physical limitations on the materials used and the 
method in which satie are constructed produce irregularities and 
rough surfaces near the leading edge of the sections. These 
imperfections, while not desirable, allow for the assumption that 


the flow is fully turbulent. This assumption precludes the 


requirement for and validation of a boundary layer transition 


model. 


II. COMPUTATIONAL APPROACH 


A. GENERAL 

Computational fluid dynamics has been largely restricted to 
well funded aerospace research and development groups due to the 
need for access to computers with sufficient power, speed and 
graphics capability. However, recent advances in computer hardware 
have eliminated the requirement for access to a super computer to 
conduct CFD studies. To demonstrate that CFD technology is mature 
enough to migrate to the field of sail design, the computational 
work for this study was conducted on Silicon Graphics Indigo and 
Iris workstations utilizing Fortran codes and Plot 3D for graphics 
presentation. This generation of workstations have largely 
redefined the conventional boundaries between workstations, 
mainframes and mini-supercomputer systems. The level of computing 
below the Unix workstation has also progressed rapidly with the 
latest generation of personal computers (PC) built around the Intel 
80486 processors. These computers probably have sufficient memory 
and speed for this application. However, at this time they are 
limited by a lack of software which has been widely available for 
the Unix operating system. It may be possible in the near future 
that new operating systems for the PC will include the necessary 
software features. The Unix based workstation, while a step above 
the PC is still well, within the resources of the sail design and 


manufacture business. 


B. SAIL SHAPE GENERATION TECHNIQUE 

To commence the study of an airfoil the first requirement 
encountered is a precise definition of the airfoil section. To 
proceed with the analysis of a sail section it is necessary to 
develop a systematic method of describing sail shapes with a 
similar high degree of precision. 

Sail design has never been an exact science. However, sail 
designers and sailors have always given attention to the same 
parameters that we see in airfoil design. The maximum camber and 
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its location aft of the leading edge are the two primary parameters 
used to describe a sail's shape (the camber of a sail is most 
commonly referred to as the “draft” by the sailing community). 
Several additional parameters are used to further describe the 
windsurfing sail section. The mast (vertical load bearing spar) 
has a cylindrical section internal to the leading edge or luff 
portion of the sail. The mast radius is directly analogous to the 
leading edge radius used in the definition of many airfoil 
sections. From the mast section aft two flat panels are used to 
fair the mast cross section to the thin segment of the section that 
compose the majority of the sail. The length of this faring has 
been defined as the luff pocket length. 

A fortran code, ‘’shape.f’ Appendix A, was developed to 
precisely define a sail section based on the geometric parameters 
mentioned above. Two models for the camber line generation have 
been included in the program. For a large number of sails the 
camber line can be approximated by a circular arc aft of the 
maximum camber location and with a second order polynomial forward. 
This model is the first camber option for sail generation which 


requires the following input arguments: 


mast radius (2.0% nominal value) 
luff pocket length 
thickness 


maximum camber or draft 
maximum camber location 
To describe a more unique section with flat segments or 
reflexed areas a higher order polynomial fit has been included as 
the second option for camber generation. For this routine the 
camber line is defined by points through which a Lagrange 
polynomial is fitted. The first three input parameters for this 
option remain the same but in place of the maximum camber and 
location the user enters n number of x,y coordinates to define the 
Camber shape. 
^Shape.f/ builds the section from the common trailing edge 


point working forward. The generation routine uses the selected 
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camber model to describe the sail camber shape aft of the luff 
dimension with an increment to add the prescribed thickness. From 
the luff dimension forward a straight line is placed tangent tHxthe 
mast radius to describe the shape. Seventy panels top and bottom 
from the trailing edge point are used to this point and a final ten 
panels are used to depict the mast radius. The 151 x,y coordinates 
that form the basic geometry of the section are written into a file 
‘shape.out. ’ 

A final segment of ’shape.f’ adds thirty wake points extended 
tangent to the trailing edge for grid generation use. The length 
between the wake points begins with the same initial dx distance 
with each subsequent length expanded by a factor of 1.1. This 
second set of data points including the wake points is written to 
‘hypgen.out’ for subsequent grid generation use. 

The primary shape used for the CFD evaluation and tunnel 
testing is depicted in Figure 2.1. This section shape was typical 
of the first fully battened sails appearing in the mid-1980’s. The 
shapes used at that time were generally highly cambered with the 
maximum camber well forward of the mid-chord point. The shape 
depicted is irregular, in particular the forward segment, when 
compared to a conventional airfoil. Obviously, a better shape from 
an aerodynamic perspective would be desirable but physical 
limitations with the mast, materials and manufacture prevent the 
building of “ideal” sections. 

Sail shapes have slowly evolved through a process of trial and 
error into shapes similar to that in Figure 2.2. The current shapes 
have less camber with the maximum camber point closer to the mid- 
cord. The shapes shown in Figures 2.1/2 were both generated using 
the circular arc and second order polynomial method. The 
“shape.out” files containing the x,y stations for the two sections 
are contained in Appendix A. 

Current interest among sail makers has been to make the aft 
sections increasingly flat. The shapes may promise to perform 
better. However, the flat areas are difficult to build into the 
sail and are prone to flutter. An advanced shape similar to the ones 
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currently being tested was derived utilizing the Lagrange 


polynomial generation method and is shown in Figure 2.3. 


C. PANEL CODE 

Panel codes represent an introductory means to model the flow 
field around a body. They are the simplest and easiest method to 
explore the flow field for the subject of this thesis. Of note, 
this family of codes can easily run solutions and be plotted on a 
modest personal computer. 

To explain the panel method we first examine Laplace's 
equation governing incompressible, irrotational, and invisid fluid 


flow: 


@ +0, =0 πο 
This expression is a homogeneous linear second order partial 
differential equation. The linearity of the equation allows the 
use of the superposition principle to describe the flow field with 
elementary flow elements. For the case of a two dimensional 
airfoil in a uniform flow the system may be expressed as a 


combination of uniform, source and vortex potential elements: 


$-0,«0,«0, ; 
where 
®_=\/_ (xcosa+ysina) 
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If the airfoil section is decomposed into a set of n panels to 
describe the section, the flow field can be represented by the 
superposition of n sources located at the mid-panel points, a 
uniform flow and a constant vortex strength. Hence, the total 
number of unknowns is n+1. 

By applying boundary conditions a solution to the system may 
be obtained. The flow tangency condition at the panel mid-points 
results in n equations. A final boundary condition is arrived at 
by applying the Kutta condition forcing the velocities of the upper 
and lower trailing edge panels to be equal. The system of 
elementary potential flows is now reduced to n+1 equations and 
unknowns which can be expressed in matrix form: 

[A] [q]=[B] 
where, A is the influence coefficient matrix, q is a column vector 
with the values of the n source strengths and the vortex strength. 
B is also a column vector equating the angular difference between 
the free stream angle of attack and the mid-point tangents and the 
wake condition. 

The coding for a panel code is fairly straightforward but 
great care must be utilized in the geometry conventions used to 
specify the influence coefficient matrix. The code ’panel.f’ 
contained in Appendix A was written by the author and adapted for 
use on sail sections. To adequately describe the irregular 
sections of interest a high number of panels was required with the 
upper limit currently 200. The airfoil shape is entered as an 
input file in one of two formats by the user. The input file can 
consist of standard x,y ordinates or be in the form of a two 
dimensional Plot-3D grid file. If a Plot-3d file is used the 
upper and lower trailing edge points along the i direction must 
also be provided. 

From a computational standpoint, the main weakness in solving 
Laplace’s equation occurs in the ’panel.f’ code when inverting the 
A matrix. The thickness of the section and number of panels 
utilized have a strong effect on the error when inverting the 
influence coefficient matrix. The actual thickness of the sail 
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section is of the order of five mil's, a very small thickness, 
which is not realistic for the code’s precision. Experience has 
shown 0.5% thickness to be a reasonable compromise. The complex 
geometries also suggest that a large number of point should be used 
to describe the section with 151 panels having been determined to 
be sufficient. 


D. BOUNDARY LAYER CODES 

The panel code provided the simplest method to obtain results 
for the pressure coefficient. Similarly, it 1s highly desirable to 
obtain viscous results in a simple and timely manner. Boundary 
layer codes represent the next logical step toward understanding 
the viscous behavior for the low Reynolds number present. 

Two boundary layer codes were explored to investigate their 
suitability for sail analysis. The first of the codes, *db12.f' a 
direct boundary layer code, was employed with very limited utility. 
The code employs a panel method routine to compute the pressure 
distribution which is then used to compute the boundary layer 
profiles. The second code utilized was the viscous invisid 
interaction code developed by Cebeci at McDonnell-Douglas Aircraft 
Company. This code carries the process used in the direct boundary 
layer code a stage further. After computing the boundary layer the 
code then adds the additional thickness of the viscous layer to the 
airfoil shape. This new effective shape is then run through the 
inviscid scheme again to compute a new boundary layer. This 
iterative routine is repeated until the change in the boundary 
layer becomes sufficiently small. The iterative boundary layer 
technique was nothing short of a total failure. This was due to 
the inability of the Smith-Hess panel routine used by the code to 
successfully handle the irregular sections. However, some insight 
into the flow separation characteristics was gained through this 
code. 

It was hoped that these two codes or methods would be of great 
value as they represent the next order of sophistication above the 


simple potential inviscid solution. Both schemes proved to be a 
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disappointment toward the analysis of the sail shapes. While these 
codes work well for standardized airfoil shapes they proved to be 
of little value for the irregular shapes. Their failure is due to 
the fact that they do not handle separated flow regions well. 
Potential flow modeling of the highly irregular sail shapes 
produced large pressure perturbations which caused the boundary 
layer to separate in the numerical solution. After several 
attempts to improve on the solutions these methods were abandoned 


in favor of the Navier-Stokes and Euler methods. 


E. EULER and NAVIER-STOKES METHOD 

The Navier-Stokes (NS) equations represent the most robust 
tool currently in use in the field of computational fluid dynamics. 
The derivation of the NS equations and their CFD solution method 
will not be discussed because they are well documented in 
References 1 and 5. The two-dimensional vector form of the 


equations may be expressed as: 
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The Euler solution is readily obtained from the NS equations when 
the viscous terms are discarded. 

The Navier-Stokes (NS) code used to examine the flow 
characteristics of the sail sections, ‘ns2.f’ was developed by 
Professor J. A. Ekaterinaris of the Navy-NASA Joint Institute of 
Aeronautics. Slight modifications to ’ns2.f’ were made in the form 
of additional write statements to save unsteady motion solutions at 
regular time intervals and to simplify steady solution inputs. A 
call to the grid rotation subroutine was added to preclude the 
input of a rotated grid when running 'ns2.f'. The major features 
of the code are: 

* Upwind differencing 

* Baldwin-Lomax turbulence modeling 

* Ability to restart the code 

Navier-Stokes solutions were obtained after approximately four 
thousand iterations from a uniform flow field, based on the density 
residuals, Figure 2.4. Solutions could be obtained with a fraction 
of the four thousand iterations if a restart was initiated from a 
previous solution (example: an 8 degree angle of attack solution is 
used to compute a ten degree solution). The NS solutions were 
computed using the Iris Indigo work stations in the CFD laboratory. 
This required approximately twelve hours of CPU time for four 
thousand iterations using a 251 by 71 grid run with a Courant 
number of 2100. While this may seem to be a huge amount of 
computer time, it must be kept in mind that this was only a small 
work station and that few of the solutions required four thousand 
iterations by using the restart feature. Euler solutions can also 
be obtained from this code in a similar manner with the correct 
switches in the input file to discard the viscous terms. 

F. GRID GENERATION 

The need for a grid system to define a flow field around a 
body arises from the necessity of transforming from the physical 
domain to a discrete computational domain. The theory of grid 
generation will not be covered in this thesis. Merkle’s text, 


Reference 5, explains this process as do several other texts in the 
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field of CFD. Grid generation for the section depicted in Figures 
2.1 and 2.2 was accomplished utilizing the code “hypgen/. 
The grid type selected for this analysis was the 'C' grid 


Navier-Stokes Convergence 
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Figure 2.4 Navier-Stokes Solution Convergence 
which is most commonly used for 2-dimensional airfoils. The ’C’ is 


in reference to the shape or manner in which this type of grid is 
wrapped around the airfoil. The i stations are along the wake and 
airfoil directions while the k stations extend from the airfoil 
surface to the far field. The primary difficulty encountered was 
ensuring the orthogonality of the grid lines around the leading 
edge and forward portion of the sail sections.’ Additional points 
along the airfoil surface were added to the luff section to smooth 
the interval distance along the i direction near the leading edge. 
For the details in the grid generation software the 'hypgen'/ user 
manual should be consulted. 

The full dimensions of the grids used are shown in Figure 2.5. 
These dimensions, 10 chords lengths ahead, above, below and 30 
chord lengths aft of the sail section, were selected to ensure that 
solutions would smoothly match the far field boundary conditions. 
The grid used for the inviscid Euler solutions is shown, Figure 
26. For the NS calculations the number of stations from the 
airfoil to the far field was increased from 41 to 71 station. This 
was done to accurately resolve the boundary layer near the airfoil 
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surface. The increased resolution can be seen in Figure 2.7. 
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III. EXPERIMENTAL APPROACH 


A. BACKGROUND 

The most frequent criticism of computational fluid dynamic 
solutions is that they may not accurately reflect the actual flow 
present. It was therefore highly desirable to validate the 
computational work conducted with experimental data for the 
sections of interest. Investigations into previous work in this 
field failed to find an adequate description of the flow around 
windsurfing sail sections. The majority of the research conducted 
in this field has been for sail and rig combinations for specific 
yacht types. The wind tunnel experiments for this thesis were 
designed with the emphasis on predicting the separation regions 


present since this is the major challenge for viscous CFD methods. 


B. WIND TUNNEL EXPERIMENTAL EQUIPMENT 

The tunnel selected for carrying out the experiment was the 32 
by 45 inch low speed wind tunnel located at NPS in Halligan Hall. 
The primary reason for the selection of this tunnel was that it is 
capable of sufficient velocities to run tests at the actual 
Reynolds numbers present in the sailing environment. The tunnel is 
also sufficiently large to allow for reasonably sized models. It 
was desired to mount the model sections vertically to take 
advantage of the large optical windows on the sides of the test 
sections for viewing tufts. In addition, by using the greater 


tunnel dimension perpendicular to the model rotation axis, the 
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blockage for the experiments was held to under twelve percent for 
a model of eighteen inch chord with up to eighteen degrees angle of 
attack. Reference 8 contains a detailed description of the tunnel 
and its use. 

The next area addressed was the choice of materials to 
construct the models from. Particular attention was given to the 
fact that the airfoils to be tested have very thin sections over 
the majority of the chord. A conventional rigid model would have 
had to be constructed of metal or fiberglass/composite materials 
for the aft section and of wood for the mast and luff portions. 
This approach was considered but abandoned due to the complexity, 
cost and long lead times in building such models. 

The other option was to assemble the model from the same 
materials as used on the actual sails. This approach simplified 
the model building process and allowed for easy modification of the 
models. The use of the same materials and construction techniques 
also precluded trying to simulate the roughness associated with the 
different fabrics, films and seams that exist on windsurfing sails. 

After deciding to use actual 2-dimensional sails as models, a 
means to hold the sail and mount the sail in the test section was 
developed. A rig was designed and built of aluminum, as is 
depicted in Figure 3.1. This rig allows for the model to be in 
tension between the top and bottom rails to maintain the correct 
shape. The leading edge or simulated mast is formed by a 0.750 inch 
diameter steel tube. This is equivalent to a mast radius of 2.08 


percent for the eighteen inch chord. The leading edge tube is 
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Figure 3.1 Wind Tunnel Model and Rig 
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freely supported and removable once the sail model is unrigged. 
This was designed to allow for the testing of several models without 
requiring a new rig. 

Mounting the rig and model in the tunnel was accomplished 
by means of a pin that runs through the upper tunnel window, 
leading edge tube and rails, fitting into a flange mounted to the 
tunnel floor, Figures 3.2 and 3.3. This method did not require 
extensive modifications to the tunnel test section. A truss type 
structure downstream of the model was incorporated to provide a 
degree of torsional stiffness to the rig. The dimension of the 
truss aft of the leading edge was arrived at to coincide with the 
breather slot downstream of the test section. This allowed for the 
linkage controlling incidence to be placed through the breather 


slots. 
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Cs WIND TUNNEL MODEL CONSTRUCTION 

The wind tunnel model that was used in the experiment was 
built and donated by Trevor Bayless at Waddell Sails in Santa Cruz. 
The materials incorporated are the same as those found on a high 
performance windsurfing sail. The leading edge luff segment was 
made of a durable dacron cloth while the aft section was made of a 
high modulus 7-mil mylar film sail material. A false seam was 
added to the leading edge to accurately pattern a seam that is 
present on actual sails. Areas beyond the body of the sail section 
were made with nylon strap material. 

To build camber into a finished sail sailmakers use a 


combination of various panel layouts and full length battens held 
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in compression to form the final shape. For the tunnel model only 
battens could be used to form the section shape since the model is 
2-dimensional. To force the airfoil shape into the desired shape 
preformed battens were used. The first set of battens, formed from 
0.125 by 0.500 inch steel, were difficult to bend into the proper 
shape and were found to be too flexible for use. This difficulty 
was overcome with the use of carbon fiber - epoxy materials. 

To build the carbon fiber battens, a female mold was first 
made of wood matching the desired camber less half the estimated 
thickness of the batten. The lamination used for the battens 


consisted of the following layers in an epoxy matrix: 


1 6 ounce S-glass fiberglass 

2 Carbon fiber reinforcing tape 
1 0.125 inch balsa wood 

2 Carbon fiber reinforcing tape 
1 6 ounce S-glass fiberglass 


The final battens proved to be very stiff and well worth the 


additional work. 


D. FLOW VISUALIZATION 

The primary means of observing the flow field around the 
section consisted of a series of tufts and a smoke stream. A 
number of different tuft materials and configurations were tried. 
Tape and cloth were placed over the lacing, used to tension the 
model, to cover the ends of the model to the rails. Flow along the 


center section was judged to be adequate for the experiments with 
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this configuration. The best tufts were determined to be very fine 
black thread with back lighting for photographs through the side 
windows. 

The smoke flow was best observed through the upper window with 
flood lighting from the side windows. Reference 8 may be consulted 
for additional information concerning flow visualization. The 
smoke wand used in the 32 by 45 inch tunnel was lengthened by 
approximately two feet to minimize the growth of the flow prior to 
reaching the test section. It was noted that the smoke flow 
behavior was sensitive to the tunnel velocity.  Velocities were 
limited to the range of 30 feet per second for this portion of the 


experiments. 
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IV. RESULTS 

A. COMPARISON OF THE COMPUTED PRESSURE DISTRIBUTIONS 

The first area examined was comparison of the computed 
pressure coefficient (Cp) distributions computed by the various 
methods. The panel method was found to be particularly sensitive 
to the thickness of the aft sections. It was desirable to keep the 
thickness value as small as possible to model the actual physical 
sections. However, the minimum feasible value was found to be half 
a percent of the chord length. The Navier-Stokes solution would 
have worked for any thickness but a value of 0.1 percent was used 
to be of the same order as in the panel method. Figure 4.1 is a 
pressure coefficient vs x/c plot from the panel code for the 11.5 


% camber section at ten degrees incidence. 


Pressure Coefficient - 10 degrees 
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Figure 4.1 Panel Method Pressure Plot 
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The abrupt changes in the Cp aft of the immediate leading edge are 
caused by the change in curvature, from the flat section at the 
luff to the second order polynomial curve at 8.5 % and at the 
transition to the circular arc section occurring at 26.5 %. The 
poor convergence of the pressure at the trailing edge is due to the 
wedge shape trailing edge. 

A comparison of the computed pressure using the panel, Euler 
and Navier-Stokes methods is depicted in Figure 4.2 for ten degrees 
incidence. The difference in the predicted Cp approximately one 
third aft of the leading edge on the lower surface using the Euler 
code is caused by the recirculating flow region. A vector plot of 
velocities for this region predicted by the Euler calculations is 
depicted in Figure 4.3 while the Navier-Stokes vector plot is shown 


in Figure 4.4. 
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Figure 4.2 Pressure Plot of Panel, NS and Euler Methods 
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Figure 4.4 Navier-Stokes Velocity Plot 
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As the incidence angle is increased the differences between 
the inviscid panel code solutions and the viscous Navier-Stokes 
predictions can be discerned. It can be seen from Figures 4.5 
though 4.7 that near stall the panel code is no longer valid. 
Figures 4.8 through 4.10 show the pressure field contours around 
the section computed by the Navier-Stokes code. The high number of 
contours around the mast radius demonstrates the importance of the 


large gradients and their effect on the stall behavior. 
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Figure 4.5 Panel and Navier-Stokes Method - 12 degrees 
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Figure 4.6 Panel and Navier-Stokes Plot - 14 degrees 
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Figure 4.7 Panel and Navier-Stokes Plot - 16 degrees 
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Figure 4.8 Navier-Stokes Pressure Contours - 12 degrees 
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Figure 4.9 Navier-Stokes Pressure Contours - 14 degrees 
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Figure 4.10 Navier-Stokes Pressure Contours - 16 degrees 


B. VELOCITY VECTOR PLOTS 

Several significant phenomena can be observed by examining the 
velocity vectors plots. At angles of attack below eight degrees a 
well defined separation and reversed flow region was predicted in 
the area near the lower luff with the Navier-Stokes code. The 
Euler code also predicted recirculating flow but estimated the 
region to be significantly larger and to be located further aft. 
The reversed flow region was first observed in the wind tunnel by 
means of the tuft survey. The size of the reversed flow region 
predicted by NS was verified by smoke flow at an angle of attack of 
four degrees. Figures 4.11 through 4.13 show this behavior by 


means of the NS prediction, smoke flow, and tufts, respectively. 
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Figure 4.13 Tuft Behavior - 4 degrees 

At angles of attack approaching stall a small separated flow 
region was predicted on the upper surface around the mast. This 
region is noteworthy because the flow quickly reattached to the 
flat luff segment and did not cause the entire upper surface to 
separate. This separation bubble was not visible during the flow 
visualization experiment, nor was this surprising. Viewing the 
bubble was highly unlikely due to its very small size. This 
prediction is pictured in Figure 4.13. 

A major goal in the flow visualization effect was to precisely 
identify the stall angle. It was observed, however, that the 
section did not stall abruptly in a manner characteristic of thin 


airfoil sections. The flow near the surface appeared to slowly 


become more turbulent until reversing around fifteen degrees. 
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Figure 4.14 Reversed Flow Region - 14 degrees 


This is an approximate angle due to the highly turbulent and random 
nature of the flow. The Navier-Stokes solutions show the boundary 
layer growing in thickness as the incidence is increased through 
fourteen degrees. This boundary layer growth is shown in Figures 
4.15, 4.17,and 4.19 with accompanying tuft photographs in Figures 


4.16, 4.18, and 4.20. 


C. UNSTEADY MOTION 

To increase the lift of sails windsurfers frequently use an 
oscillatory motion to increase the maximum lift generated by the 
Sail. An effort to understand this phenomena was made using the 
unsteady options in 'ns2.f'. Several simulations were made with 


the test airfoil but only one solution will be presented. To simulate 
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Figure 4.18 Upper surface - 14 degrees 
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the sailors actions a sinusoidal motion of four degree about a mean 
incidence of fourteen degrees was used. The reduced frequency was 
estimated to be 0.19 for an oscillation with a 1.5 second periods 
The code is started from a steady solution for time equal zero with 
a frequency shift of minus one half pi. The motion simulated is 
shown in Figure 4.21 as a function of time. The velocity field 
that results from the motion is shown in Figure 4.22 through 4.26. 

The velocity vectors for the simulated motion show attached 
flow beyond the steady motion stall angle. The flow in fact stays 
attached through eighteen degrees but separates after the motion 


starts down. 
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V. CONCLUSIONS 


In the computational and experimental investigation of the 
flow field about windsurfing sail sections several important 
findings were made. The first area of surprise was the extent of 
flow separation that is present over several areas on the sections. 
This was evident in both the Navier-Stokes calculations and in the 
flow visualization experiments. Moreover, the degree of separation 
was the primary reason for the failure of the two boundary layer 
methods that were applied to this problem. While the panel method 
represents the most elementary code used on the problem it compared 
favorably with the more costly inviscid Euler routine. Since the 
Euler code is nearly as expensive in computer time as the Navier- 
Stokes method its use can not be justified. Instead, much useful 
information can be obtained from the Navier-Stokes solutions. 

The Navier-Stokes code, ’ns2.f’, while performing well, has 
been improved upon. Specifically, a variable time stepping routine 
has been added which offers a reduction in the number of iterations 
and computer time by over one half. Only one turbulence model, 
Baldwin-Lomax, was incorporated in the version that was evaluated 
in the thesis. Several other models have since been added and 
could be evaluated. Further investigations of the effect of sail 
pumping used in competitive board sailing would be of high interest 
to sailors and to aerodynamicists. Unfortunately, a failure of the 
computer system in the CFD laboratory prevented the modeling of a 


realistic ramp motion for use in 'pumping' simulations. 
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The experiments undertaken achieved their primary objective to 
locate the separation regions. However, the visualization 
experiments could be greatly improved upon. This could be 
accomplished with the implementation of the laser sheet technique. 
The rig and model could also be modified to include the examination 
of unsteady motion effects. 

Finally, it is clear that advanced CFD methods, such as 
Navier-Stokes solvers, can be used in this field. Sail design is 
a field that thrives on new technology. The advances in computer 
technology make it only a matter of time when CFD techniques will 


migrate to the field of sail design. 
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APPENDIX A: COMPUTER PROGRAMS 
Some of the computer programs that were used during the 
research are presented in the following section. This listing is 
by no means a complete listing of the software utilized. Not 
listed but of significant utility are 'hypgen' and 'plot3d' both of 


which are well documented. 
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*****  shape.f 

* Lt Matthew Avila 

* Thesis Advisor: Prof Platzer 

* Sail section generation subroutine 

* Parabolic/circular arc or Lagrange Interpolation 

* to define camber line 

* Wake points added for hypgen use ( hypgen.dat ) 
* X-Y data points ( shape.out ) 

* 


DIMENSION X(0:150),Y(0:150),XNODE (30) , ΥΝΟΡΕ (30) 
DIMENSION xx(300),yy(300) ,zz(300) ,Xwake (30) , Ywake (30) 
OPEN (101, FILE=’ shape.out’ , STATUS=’ UNKNOWN’ ) 

OPEN (102, FILE=’hypgen.dat’ , FORM=’ formatted’ ) 
pi=3.14159265 

X(0)=1.0 

Y(0)=0.0 

X(150)=1.0 

Y(150)=0.0 


PRINT *,'This program generates X-Y section coordinates and' 

PRINT *,'Plot3d file for with the HYPGEN grid generation 
program.’ 

PRINT *,'151 points are created for the X-Y plot and a 
2 

PRINT *,’grid in plot3d format.’ 


10 PRINT *,’Enter camber line generation method’ 

PRINT *,’ 1 = parabola forward section / circular arc 
aft’ 

PRINT *,’ 2 = Lagrange polynomial defined by user’ 


READ *, TYPE 


IF (TYPE.EQ.1) THEN 


GOTO 20 
ELSEIF (TYPE.EQ.2) THEN 
GOTO 50 
ELSE 
PRINT *,'You must enter 1l or 2' 
GOTO 10 
ENDIF 
* Input Variables for Parabolic-arc method 
20 PRINT *,'Enter mast radius (.02 normal)' 
READ *,r 


PRINT *,’Enter Maximum camber’ 

READ *, Ymax 

PRINT *,’Enter Maximum camber location’ 
READ *, Xmax 

PRINT *,‘Enter Luff Pocket Length’ 

READ *,Xluff 

PRINT *,’Enter thickness’ 

READ *,t 


50 


tzt/2 

nfwdsINT (Xluff*70) 

naftz70-nfwd 
rl-((1-Xmax)**24«Ymax**2) /2/Ymax 


* Aft sections / Parabolic-Arc 
dxaft=(1-Xluff)/naft 
DO 30 I=1,naft 
Xl=I*dxaft 
IF (X1.LE.(1-Xmax)) THEN 
SisASIN((1-Xmax-X1)/r1) 
Yl=r1*COS (Si) +Ymax-r1 
ELSE 
Yl=Ymax- Ymax* ( (X1+Xmax-1) **2) / (Xmax**2) 
ENDIF 
X(1)=1.0-X1 
Y(I)sYl-t 
X(150-I)z1.0-X1 
Y(150-I)=Yl+t 
30 CONTINUE 


* Wake section / parabolic-Arc 
Si = ASIN( (i.-Xmax) / (ri-Ymax) ) 
Xwake (1)=1+dxaft 
Ywake (1) =-dxaft*TAN (Si) 
DO 40 I=2,30 
dxaft=dxaft*1.2 
Xwake (I) =Xwake (I-1)+dxaft 


Ywake(I)= - ( Xwake(I) -1. )*TAN(Si) 
40 CONTINUE 
GOTO 100 
* Input Variables for the Lagrange method 
50 PRINT *,’Enter number of nodes to describe camber line’ 
READ *,nn 


DO 60 I=1,nn 
PRINT *, “Enter (x,y) for node', 1I 
READ *,XNODE(I) , YNODE (TI) 
60 CONTINUE 


PRINT *,’Enter mast radius (.02 normal)’ 
READ *,r 

PRINT *,’Enter Luff Pocket Length’ 

READ *,Xluff 

PRINT *,'Enter thickness' 

READ *,t 

t=t/2 

nfwdsINT (Xluff*70) 

naft=70-nfwd 


* Aft sections / Lagrange method 
dxaft=(1-Xluff) /naft 
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DO 70 1-1 naft 
دع لاس لګ‎ 1171٢ 
Y1=0 
DO 8O0 ٢٣٢ pn 
Z = YNODE (K) 
DO 90 ΤΕΙ, ΠΗ 
IF (K.NE.L) Z=Z*(X1-XNODE (L) ) / (XNODE (K) -ΧΝΟΡΕ (1) ) 


90 CONTINUE 
Y1 = Y1 + Z 

80 CONTINUE 

X(T) =I 

Y (1) = YI L 

x(150-I) = X1 

Y(150-I) = Yl+t 
10 CONTINUE 


* Wake section / Lagrange 
100 Si = ATAN( (Y (1) +Y(149))/2/dxaft) 
Xwake (1) =1+dxaft 
Ywake (1) =-dxaft*TAN (Si) 
DO 110 I=2,30 
dxaft=dxaft*1.2 
Xwake (1) =Xwake (1-1) +dxaft 
Ywake(I)= - ( Xwake(I) -1. )*TAN(Si) 
110 CONTINUE 


* Lower Luff Section / Common for both methods 
Thel-pi/2-ATAN(Y(naft)/(Xluff-r)) 
X (70) =r+r*COS (thel) 
Y (70) =-r*SIN(Thel) 
dx=X (70) -X(naft) 
slope-(Y(70)-Y(naft))/dx 
dxfwdsdx/nfwd 
DO 120 I-21, (nfwd-1) 
X (I+naft)=Xluff+dxfwd*1I 
Y (I+naft)=Y (naft)+(X(I+naft) -X(naft))*slope 
120 CONTINUE 


* Upper Luff Section / Common both methods 
The2=pi/2-ATAN( (Xluff-r) /Y(150-naft) ) 
X (80) =r-r*SIN (The2) 
Y (80) =r*COS (The2) 
dx=X(150-naft) -X(80) 
slope=(Y(150-naft) -Y(80) ) /dx 
dxfwd=dx/nfwd 
DO 130 Iz1, (nfwd-1) 
X (1+80) =X (80) «dxfwd*I 
Y (1+80) =Y (80) + (X(1+80) -X (80) ) *slope 
130 CONTINUE 
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* Leading Edge Section 
Theda=1.50*pi-Thel-The2 
dTheda=Theda/10 
DO 140 I=1,9 

X(70+I)=r+r*CoS (dTheda*I+Thel) 
Y (70+1) =-r*SIN(dTheda*I+Thel) 
140 CONTINUE 


* SHAPE.OUT ROUTINE 
DO 150 I=0,150 
WRITE (101,*)X(1),Y (1) 
150 CONTINUE 
* HYPGEN ROUTINE 
211 


1 
1 


im 
jm 
km 


l H | 


j=0 

do 160 i=30,1,-1 

j=j+1 

xx(j)=Xwake (i) 

yy (3) =0. 

zz(j)sYwake(i) 
160 continue 


do 170 is0,150 

ii = i431 

xx(ii) = x(i) 

VV (11) =-O. 

ZZ (11) = y(1) 
170 continue 


do 180 i=1,30 

xx (1+181)=Xwake (i) 

yy (i+181)=0. 

zz(i+181)=Ywake (i) 
180 continue 


rewind 102 

write (102,*) 1 

write (102,*) im,jm,km 

write (102,*) (xx(i), i=1,211), 
> (yy (i), i=1,211), 
> (zz(i), i=1,211) 
close (102) 


END 
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*****  Danel.f 
* Lt Matthew Avila Last modified 19 Oct 92 
* Thesis Advisor: Prof Platzer 
* Panel Method Program 
* Airfoil shape input from Plot3d 2 dimensional grid or 
* X - Y data file 
* پ پ ې‎ xy 
PARAMETER (ID=300,JD=80) 
COMMON /SET1/X(ID) ,Y(ID) ,DIST(ID) , SINT (ID) , COST (ID) 
COMMON /SET2/A(ID, ID) ,P1I,NPANEL, NODES, NTOT 
DIMENSION XG(ID,JD) , YG(ID,JD) 
CHARACTER*80 FILEIN 
INTEGER TYPE 
1 FORMAT (A) 
PI=3.1415926585 


10 PRINT *, "ENTER DATA TYPE’ 
ΕΕ 1 PLOT3D FILE’ 
PRINT *, ' 2 AX c Y DATA: 
READ (*,*)TYPE 
IF (TYPE.EQ.1) THEN 

GOTO 20 
ELSEIF (TYPE.EQ.2) THEN 
GOTO 40 
ELSE 
PRINT *,'YOU MUST ENTER 1 OR 2' 
GOTO 10 
ENDIF 


* READ PLOT3D GRID POINTS FOR AIRFOIL 
20 PRINT *, ’ENTER 2-D GRID FILE NAME’ 
READ (*,1) FILEIN 
OPEN (UNIT=201, FILE=FILEIN, STATUS=’ OLD’ , FORM=’ FORMATTED’ ) 
PRINT *,’ENTER I FOR TRAILING EDGE LOWER, UPPER’ 
READ *,ITEL,ITEU 
REWIND 201 
READ (201,*)IMAX,JMAX 
PRINT *,IMAX,JMAX 
READ (201,*) ((XG(I,J) ,I=1, IMAX) ,J=1,JMAX), 


> ((YG(I,J),I=1,IMAX), J=1,JMAX) 
I = O 

DO 30 JsITEL,ITEU 

I=1+]1 1 


X (I) =XG(J,1) 
Y (I) =YG(J,1) 
30 CONTINUE 


NPANEL = ITEU - ITEL 
NODES = NPANEL + 1 
NTOT = NODES + 1 
GOTO 60 
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* READ X-Y DATA POINTS FROM A FILE 
40 PRINT *,'ENTER X-Y FILE FOR AIRFOIL' 
READ (*,1)FILEIN 
OPEN (UNIT=201, FILE=FILEIN, STATUS=' OLD’ , FORM=’ FORMATTED’ ) 
PRINT *, ’ENTER NUMBER OF DATA POINTS’ 
READ *,NODES 
NPANEL NODES - 1 
NTOT NODES + 1 


DO 50 I=1,NODES 
READ (201,*)X(1),Y (1) 
50 CONTINUE 


* ANGLE OF ATTACK FOR COMPUTATION ENTERED 

60 PRINT *,'Enter Angle of attack (degrees)' 
READ *,AOA 
SINA=SIN(AOA*PI/180) 
COSAZCOS (AOA*PI/180) 


* SLOPE AND LENGTH OF PANELS COMPUTED 
DO 70 J=1,NPANEL 
DX=X (J+1) -X(J) 
DY=Y (J+1) -Y (Jd) 
DIST (J) =SQRT (DX**24+DY**2) 
SINT (J) =DY/DIST (J) 
COST (J) =DX/DIST (J) 
70 CONTINUE 


* EXECUTE INFLUENCE COEFFICIENT, GAUSS AND VELOCITY SUBROUTINES 
PRINT *,’PANELS COMPUTED’ 
CALL COEFF (SINA, COSA) 
PRINT *,’COEFFICIENTS MATRIX COMPUTED’ 
CALL GAUSS (SINA) 
PRINT *,’GAUSS REDUCTION COMPLETE’ 
CALL VELOCITY (SINA, COSA) 
PRINT *,’ PRESSURE COEFF. COMPUTED (panel.out) ’ 
END 
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INFLUENCE COEFFICIENT SUBROUTINE‏ ٭ 
SUBROUTINE COEFF(SINA,COSA)‏ 
PARAMETER (ID=300)‏ 
COMMON /SET1/X(ID),Y(ID),DIST(ID),SINT(ID), COST (ID)‏ 
COMMON /SET2/A(ID,ID),PI,NPANEL,NODES,NTOT‏ 


DO 110 K=1,NODES 
A(NODES,K)=0.0 
110 CONTINUE 


* MIDPOINT NORMAL VELOCITY=0 
* ['th PANEL MIDPOINT 
DO 120 I=1,NPANEL 
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XMID=0.5* (X (1) +X(1+1)) 
YMID=0.5* (Y (1) +Y (1+1)) 
A(I,NODES)=0.0 


* J 'th PANEL 
DO 130 J=1,NPANEL 

DLN = 0.0 

BETA = PI 
* I AND J THE SAME 

IF (J.E0.1) GO TO 140 

DX1 = XMID-X(J) 

DX2 = XMID-X(J+1) 

DY1 = YMID-Y(J) 

DY2 = YMID-Y(J+1) 

D1 - SQRT(DX1**24DY1**2) 

D2 - SQRT(DX2**24DY2**2) 

DLN = LOG(D2/D1) 

BETA=ATAN2 (DY2*DX1-DX2*DY1,DX2*DX1+DY2*DY1) 
140 SINTIJ-SINT(I)*COST(J)-COST(I)*SINT(J) 


COSTIJ-COST(I)*COST(J) *SINT(I)*SINT(J) 
A(I,J)=(BETA*COSTIJ+DLN*SINTIJ)/PI/2 
B- (DLN*COSTIJ-BETA*SINTIJ)/PI/2 
A(I,NODES) =A(I,NODES) +B 
* KUTTA CONDITION APPLIED TO FIRST/LAST PANEL 
IF ((I.GT.1). AND. (1.LT.NPANEL) ۷٤٢ 
A(NODES,J)=A(NODES,J) -B 
A (NODES , NODES) =A(NODES, NODES) +A(I,J) 
130 CONTINUE 
A(I,NTOT) =SINT (I) *COSA-COST (I) #SINA 
120 CONTINUE 
TEMP2- (COST (1) «COST (NPANEL) ) *COSA- (SINT (1) «SINT (NPANEL) ) *SINA 
A (NODES , NTOT) TEMP 
END 
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* GAUSS SUBROUTINE TO SOLVE COEFFICIENT MATRIX 
SUBROUTINE GAUSS (*) 
PARAMETER (ID=300) 
COMMON /SET1/X(ID),Y(ID),DIST(ID),SINT(ID), COST(ID) 
COMMON /SET2/A(ID,ID),PI,NPANEL,NODES, NTOT 


* SEARCH FOR THE LARGEST PIVOT 
DO 200 I=2,NTOT 
IMsI-1 
IMAX= IM 
AMAX=ABS (A (IM, IM) ) 
DO 210 J=JI,NODES 
IF (AMAX.GE.ABS(A(J,IM))) GO TO 210 
IMAX=J 
AMAX-ABS (A(J,IM)) 
210 CONTINUE 
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* SWITCH ROWS IF NECESSARY 
IF (IMAX.EQ.IM) GO TO 220 


DO 230 J = IM, NTOT 
A(IM,J) = A(IMAX,J) 


A (IMAX,J) =TEMP 
230 CONTINUE 


* REDUCE BELOW THE PIVOT DIAGONAL 
220 DO 250 J=I,NODES 


R = A(J,IM) /A(IM, IM) 
DO 240 K=IM,NTOT 
A(J,K) = A(J,K)-R*A(IM,K) 
240 CONTINUE 


250 CONTINUE 


* REDUCE ABOVE THE PIVOT DIAGONAL 
500۰۰260٣۰0  - 77 


R = A(J,IM)/A(IM, IM) 
DO 270 K=IM,NTOT 
A(J,K) = A(J,K)-R*A(IM, Κ) 
270 CONTINUE 


260 CONTINUE 
200 CONTINUE 


* CHANGE PIVOT VALUES TO 1.0 
DO 280 I=1,NODES 
A(I,NTOT)-A(I,NTOT)/A(I,I) 
A(I,I)=1.0 

280 CONTINUE 
END 
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* VELOCITY DISTRIBUTION SUBROUTINE 
SUBROUTINE VELOCITY (SINA, COSA) 
PARAMETER (ID=300) 
COMMON /SET1/X(ID),Y(ID),DIST(ID),SINT(ID),COST(ID) 
COMMON /SET2/A(ID,ID),PI,NPANEL,NODES, NTOT 
DIMENSION CP(ID),Q(ID) 
OPEN (301,FILE-'panel.out',STATUS-'unknown') 
YMULT=20.0 


DO 310 I=1,NPANEL 
310 Q(I)=A(I,NTOT) 
GAMMA=A (NODES,NTOT) 


* I'th PANEL MIDPOINT 
DO 330 I=1,NPANEL 
XMID=0.5*(X(1)+X(1+1)) 
YMID=0.5*(Y(I)+Y(1+1)) 
VTAN=COSA#COST (I) +SINA*SINT (I) 
* J'th PANEL 
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300 


320 


330 


CON 
END 


DO 320 J=1,NPANEL 


DLN=0.0 

BETA=PI 

IF (J.EQ.I) GO TO 300 

DX1 = XMID-X(J) 

DX2 = XMID-X(J+1) 

DY1 = YMID-Y(J) 

DY2 = YMID-Y(J+1) 

D1 = SORT (DX1**24+DY1**2) 
D2 = SORT (DX2**2+DY2**2) 
DLN = LOG(D2/D1) 


BETA=ATAN2 (DY2*DX1-DX2*DY1,DX2*DX1+DY2*DY1) 
SINTIJ=SINT(1)*COST (J) -COST(1) *SINT (J) 
COSTIJ=COST (I) *COST(J) +SINT(I)*SINT(J) 

AA- (BETA*COSTIJ«DLN*SINTIJ)/PI/2 

B» (DLN*COSTIJ-BETA*SINTIJ)/PI/2 
VTAN-VTAN-B*Q(J) «GAMMA*AA 

CONTINUE 

CP(1)=1.0-VTAN**2 

WRITE (301,*)XMID,CP(I) 

TINUE 
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******  rotate.f 

* Lt Matthew Avila 

* Thesis Advisor: Prof Platzer 

* Plot3d grid Rotation Program 

Aki kk k 
parameter (1d=300,kd=80) 
dimension X(id,kd),Z(id,kd) 
character*15 filein,fileout 


* Read in original grid 
print *,'Enter plot3d grid file' 
10 format (A) 
read (*,10)filein 
open (20,file=filein, form=’ formatted’ ) 
rewind 20 
read (20,*) imax,kmax 
read (20,*) ((X(i,k),1=1, imax) ,k=1,kmax), 
> ((Z(i,k),i1s1,imax),ks1,kmax) 


* Rotate the grid points 
print *,'Enter the angle to rotate the grid by' 
read (*,*)theta 
theta=theta*3.1415/180 
cost=cos (theta) 
Sint=-sin(theta) 
do 30 i-1,imax 
do 30 k=1,kmax 
xold=X (i,k) 
zold=Z(1,K) 
X(i,k)sxold*cost - zold*sint 
Z(1,k)=zold*cost + xold*sint 
30 continue 
* Write output file 
print *,'Enter filename for rotated grid’ 
read (*,10)fileout 
open (40,file-fileout,form-'formatted') 
rewind 40 
write (40,*) imax,kmax 
write (40,*) (( X(1,k), i=1,imax), k=1,kmax), 
> (( 2(1,k), i=1,imax), k=1,kmax) 


stop 
end 
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***** ns2.in - input file for ns2.f 


MACH,  ALFAO, ALFAl, ALFARE, REDFRE, REYNOLDS 


0.10, 14.00, 4.0, O00, 0.19, 0.80 
ED2X, ED2Y, ED4X, ED4Y, ED 
0.00, 0.00, 0.030, 0.030, 220 
DT, COUR, NITER, NEWTIT 
0.0002, 2100., 2160, 1 
RSTRT,  OSCIL, RAMP, NPER TSHIFT 
true, true, false, 1440, ۔‎ 5 
TIMEAC, IMPLBC, EXPLBC  CIRCOR 
1, false, true, false 
VISC, BLTM, JKTM, RNGTM 
true, true, false, false 
ITEL, ITEU 
31. 221 
UNSTST, NFRAME 
true, 8, 
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Mach Í Free stream Mach number 
Alfao0 Σ Agle of attack, also mean angel of attack for unsteady 
Alfal : Amplitude of Oscillatory motion 
Redfre : Reduced frequency k - omege * c / 2U 
Reynolds : Reynolds number Re = cU/n 
ED2x Ξ X-direction 2nd order explicit smoothing ( e2x = 0.00 
subsonic, 

0.25 < e2x < 0.50 
transonic ) 
ED2z : Z-direction 2nd order explicit smoothing ( e2z - 0.00 
subsonic, 

0.10 < e2z < 0.20 
transonic ) 
ED4x : X-direction 4nd order explicit smoothing ( e4x = 0.03 
subsonic, 

e4x = 0.05 
transonic ) 
ED4z : Z-direction 4nd order explicit smoothing ( e4z = 0.03 
subsonic, 
e4z 5 

transonic ) 
ED : Scaling of Implicit smoothing 
ISPEC : Spectral radious parameter 
Dt : Time step 
Cour : Courant Number Cu - dt * L max 
Niter : Number of Iteration in this run 
Newtit : Newton subiteration within each timestep 
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RSTRT : Restart 


OSCIL : Oscillatory motion A(t) = AO + Al * sin ( Kk *k M* t) 
RAMP : Ramp motion 

NPER : Number of time steps in one period of oscillation, 
dt-T/Nper 

TSHIFT : Time shift in radiats to start oscilation for any a(t) 
TIMEAC : Time accureat Tacc=1 and for Jacobian Scaled Dt, 
Tacc=0 

IMPLBC : Implicit wall bc Treatment 

EXPLBC : Explicit wall bc Treatment 

VISC : Viscid or inviscid Boundary Cconditions (if false 
inviscid) 

TURBL : Boldwin-Lomax eddy viscosity (if false laminar only) 
JKTM : Johnson-King eddy viscosity (not available) 

RNGTM : RNG eddy viscosity (not available) 


ITEL, ITEU : Lower and Upper trailing edge I locations 


UNSTST (Set true for use with unsteady motion starting from an 
steady 


state. This will initialize time = O for the unsteady 
motion 

that starts from a steady state restart. 
NFRAME : Number of solutions saved for one cycle of 
oscillation. 


Saved starting at unit 41. 


Read grid from unit fort.11 and the flow from unit fort.31 
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